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Mettio4 and Apparatus for Qua^ 

The present invention relates to a method and apparatus for quantifying material 
properties. More particularly, the invention relates to such methods^and apparatus 
which rely upon a spectral analysis of light remitted, emitted, and/or transmitted from 
the rnaterial or object under test. . 

There are many situations in which non-invasive arid- non-destructive methods of 
analysing an object or material are desirable or even necessary. For example, in the 
medical field, it is. often desirable to be able to analyse the properties of tissue in vivo, 
avoiding the need for biopsies and other invasive and destructive techniques. One 
technique employed for this purpose relies upon an analysis of light reflected (or more 
correctly "remitted") by an area of tissue. Assmning that the illuminating Hght 
penetrates the tissue to sdme depth and is reflected (or scattered or/and absorbed) to 
different degrees at different depths, ari analysis of the remitted light can provide useful 
information about the properties of the tissue. 

Techniques such as that described in the preceding paragraph build upon so-called 
"physics-based" approaches to computer vision. In contrast to traditional image 
analysis techniques which analyse a final image, physics-based approaches rely on\ 
fundamental understanding of the image formation process. This process consists of an 
involved series of interactions of the iUuminating Ught at. and wifliin. object surfaces. 
By using an understanding of these interactions, relationships between scene parameters 
and image values are estabUshed. Physics-based approaches have been shown to be 
successful in many different areas of image understanding, such as obtaining shape 
from shading, and removing the effects of variations in illumination prior to image 
segmentation. 

A technique for non-invasively analysing skin stiiicture is described in International 
patent appUcation pubUcation number WO98/22023. This techniques uses a physics- 
based computer vision approach. Tiie technique described in this publication makes use ' 
of a caUbration step in order to take into account flie fact tiiat die skin iUuinination 
conditions may vary from one measurement to the next. The caUbfation step, j^uires. 



the use of some additional equipment and also adds to the overall- analysis time. Whilst 
this may not present a problem for analysing sMn, it may present a problem for 
analysing other body tissues. For example, it may not be possible to provide additional 
calibration components on probes used, for" internal examination, (such as endoscopes), 
whilst calibration is difficult to achieve when analysing the fundus of the eye due to. the 
complexity of the eye stmcture and the unpredictabiUty.of the intervening lens. It will 
be appreciated that calibration may also present problems when analysing non-body 
materials where proximity or contact of the calibration medium with the material 
. analysed is difficult or impossible. 

Substantial work has been carried out to develop image analysis algorithms which are 
able to identify different objects irrespective of the illuminating li^t. Many of the 
techniques developed are based around the linear model of surface reflectance as 
proposed in L; Maloney and B. Wandell, "Color constancy: a method for recovering 
surface spectral reflectance", J. Opt. Soc. Am. A 3, 29-33 (1986). This approach is 
based on the idea that the surface reflectance of any object within an imaged scene can 
be expressed as a weighted sum of basis spectral reflectance functions: 

and that the illuminating U^t can sunilarly be expressed as a weighted sum of basis 
lights. It has been shown that only a small number of basis functions are required to 
obtain accurate approximations to the surface reflectances of many naturally occurring 
objects and also the spectral variation of natural dayli^t. 

With this technique it is possible to recover the vector of weighting constants ctj from a 
vector of image values and thus specify the spectral reflectance of the imaged object at 
every pixel. Every potential imaged object/object characteristio will have a unique 
spectral reflectance. Thus, if the spectral reflectance can be determined using a Unear 
model, then the parameter vector can be specified. With this approach it should be 
possible to recover a parameter vector from the vector of image values at each pixel. 
Unfortunately the method is only able to recover the weighting constants Oj to within a 
multipUcative scaling factor and thus cannot be used to specify the exact spectral 
reflectance and therefore ttie exact parameter vector. 



An approach to geometry-insensitive segmentation of images has been developed in G. 
Healey, 'TJsing colour for geometry-insensitive segmentation," J. OptrSoc. Am. A 6, 
920-937 (1989), and is based on the idea, of noraialised colour. With Ihis approach 
image values are first divided by an estimate of normalised colour. This estimate is 
based on approximating the incoroing signal from colour pixel values by using a finite- 
dimensional linear approximation to represent the colour signal. Using these 
normalised values, different metal and dielectric materials can be identified across an 
imaged scene in which the geometry varies considerably. 

A similar technique has been applied to evaluate bum injuries {M. A. Afiromowitz, G. 
S. van liew and D. M. Heimbach, "Clinical evaluation of bum injuries using an optical 
reflectance technique", IEEE trans. Biomed. Eng. BME-34, 114-127 (1987), and M. A. 
Afromowitz, J. B. Callis D. M. Heimbach, L. A. Desoto and M. K. Norton, 
"Mulitspectral imaging of bum v^Qunds: a new clinical instrument for evaluating bum 
depth", IEEE tran. Biomed. Eng. 35, 842-849 (1988)}. In this case, RGB image values 
were normalised by dividing them by the response of an IR filter. From the nomialised 
values it was possible to assess the extent of bum damage across a given area of imaged 
skin. 

There exists a need for a non-invasive technique for analysing an object or material 
(which may be complex, for example multi-component and/or multi-layer and which 
may be solid, gaseous, liquid, etc) which does not require calibration to take into 
account changing illumination conditions. It is an object of the present invention to 
provide such a technique. This and other objects are achieved by designing optimal 
spectral filters which will select light components which can be mapped uniquely to the 
object/material properties of interest, with a sufficiently high resolution. 

According to a first aspect of the present invention there is provided a method of 
deteniiiriing a property or properties of each of a set of filters, which filters are used to 
select specific wavelength ranges in a system which relies upon a spectral analysis of 
remitted, emitted, and/or transmitted light to quantify a parameter or parameters of an 
object or naaterial, the method comprising the steps of: 
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1) defining a set of potential filter properties 

.2) defining one or more image quotients, the or each image quotient for a region 
being, obtained by dividing the quantified output of a given filter for "that region, the 
"image value" for that filter, by another image value for that same region; " 

3) for an object or material to be analysed and for said defined set of potential 
filter properties and for said given image quotients, obtaining a function mapping points 
in parameter space to points in image quotient space; 

4) deteimining whether the mapping function provides a 1:1 correspondence 
between points in parameter space and points in image quotient space; and 

5) if the mapping function does not provide a 1:1 correspondence, rejecting said 
potential filter properties, repeating steps 1) to 4) and, if the mapping function does 
provide a 1:1 correspondence accepting the potential filter properties as a candidate set 
of filter properties. 

The key step in the present invention is that of identifying a set of filter properties and 
image quotients and a mapping function which maps image quotients to object/material 
parameters with , a 1:1 correspondence. This first step may require consideration of 
many potential filter properties and corresponding mapping functions. 

Preferably, for each filter the method of the present invention is used to determine the 
centre wavelength of the filter. The method may additionally be used to determine the 
full width half maximum (FWHM) of the filter characteristic. 

Preferably, step 3) comprises constructing a Jacobian matrix for the mapping function 
with respect to said parameter(s), and obtairdng the determinant of that nciatrix. If the 
determinant is strictly positive or stricfly negative over the entire parameter space, then 
a 1:1 correspondence between points in parameter space and points in image quotient 
space is assumed. It will be appreciated that these operations may be performed using 
multiple processing steps or may be combined into a single processing step. 

Embodiments of the present invention may involve the calculation of an error, 
indicative of the accuracy of parameter recovery obtained using said mapping function. 
The error may be calculated as follows: 



,* a) calculate the error associated witti image acquisition for each vector ; 
of each image quotient; 

b) from the iniage quotient vector error, calculate the maximum possible error 
in each component of the parameter vector across the whole of parameter space; 
and 

c) use the vector of parameter errors at each point within parameter space to . 
measure the accuracy of parameter recovery. 

The present invention may be incorporated into many different filter property 
calculation schemes. For example, in a scheme using a genetic algorithm, the method 
may be used to identify a pluraUty of candidate filter parameter sets. The method of the 
present invention is then applied repeatedly to find an optimal filter parameter set using 
the candidates. Alternatively, the -method may be employed in a scheme using a 
gradient descent algorithm. In such a scheme, the method of the present invention is 
employed to identify a first candidate set of filter parameters. This set is then used to 
make a selection of another set of filter properties, and the process repeated as often as 
necessary to arrive at an optimal solution. Of course, any suitable optimization 
algorithm can be used to compute an optimal solution or a solution which has sufficient 
accuracy. 

The present invention may be employed in many different technical fields. However, 
one field in which it is particularly ajpplicable is that of medical diagnosis and testing, 
for example analysing skin structure. 

Each time the method steps 1) to 4) are carried out, the image quotients may or may not 
be changed. That is to say that, for each repetition, the potential filter properties and the 
image quotients may be changed, or only the potential filter properties- may be changed. 

According to a second aspect of the present invention there is provided apparatus for 
analysing an object or material having means for conducting a spectral analysis of 
remitted, emitted, and/or transmitted light to quantify a parameter or parameters of an 
object or material, the apparatus comprising a plurality of filters for sphtting said hght 
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into respective components, the filters having properties obtained by using the method 
of the above first aspect of the invention. 

It wiU be appreciated that the filters of the apparatus may be implemented in optical, 
electrical/electronic, or software form. 

For a better understanding of the present invention and in order to show how tiie same 
may be carried into effect, reference wUl now be made by way of example to the 
accompanying drawings, in which: 

Figure 1 illustrates a function f which maps an image quotient i to a material parameter 

p; 

Figure 2 illustrates a combination of. errors in parameter spape, mapped to image 
quotient vector space; 

Figure 3 illustrates a model of the layered structure of skin; 

Figures 4a and 4b show how the reemitted spectrum (intensity vs wavelength) varies for 
different melanin and blood levels respectively; 

Figure 5 illustrates a set of filters suitable for analysing blood and melanin levels in 
skin; and 

Figure 6 is a flow diagram illustrating a method of defining a set of filter properties for 
use in analysing the properties of an object or material. 

Most physics-based approaches used in computer vision employ what is known as the 
dichromatic model. A detailed explanation of this model is given in S.A. Shafer, 
"Using colour to separate reflection components," Col. Res. Appl. 10, 210-218 (1985). 
The dichromatic model states that light remitted from an object is the sum of two 
components, a "body" component and a "surface" component. The body component 
relates to physical processes occurring after penetration of light into the material and the 
surface component relates to reflections which take place at the surface of the object. 
The body component is a function of the spectral characteristics of tiie object, whereas 
• the sijrface component depends only on the object geometry and the incident light. The 
model states further that each component can be considered as die product of a 
geometrical tann and a wavelength dependent term. 



Jh a typical analysis system, light remitted from an object can be measured using a 
digital camera with a .small number of optical filters. Image values - brightness or 
"intensity" - for each image location (x, y) for a"given filter (the nth filter) are given by: 

^ KtjB(X)S{X)R'{X)dX + KsjB{X)Ii^(X^ (2) 

where K\y and Ks are the geometric terms of the body and surface components 
respectively and and C" are colour terms. The first integral in equation (2) is the 
product of three terms: E(A) is the illuminating light, S(Ji) is the spectral remittance 
from the body of the imaged object, and is the spectral response of the nth optical 
filter. In the second integral there are. only two terms as there is no wavelength 
dependence on the surface component of reflection. The dichromatic reflection model 
is very important for 3-D scene analysis as it allows for both colour and geometrical 
analysis of objects within a scene. 

The approach presented here is concerned with problems in which the surface term in 
equation (2) can be neglected. Although this would appear to limit the class of 
problems to scenes containing objects which do not have a specular component of 
reflection, such as many organic objects, image processing algorithms have been 
developed to allow removal of this component of reflection giving greater applicability 
of the technique. For example, it has been shown that it is possible to remove the 
highlights from complex images containing inhomogeneous dielectrics. It is also 
possible to" remove the surface component of reflection using a polarising filter. Once 
this component has been removed from image* data, it will be possible to use the 
techniques described here. 

A key issue is to show that the technique proposed here is valid for problems where the 
intensity of the illuminating light is imknown (whilst assuming that the spectral 
definition of the illuminating light is known). For this purpose the incident light is 
written as: 

B(X)^eoEo(X) (3) 
where e© is a wavelength independent scaling factor. Equation (2) now becomes 
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t^{x,y) = €jEoiX)S{X)II^{X)dX * (4) 

where e = eoKb.. A digital camera records an iVrdimensional vector of image values at 
-each location (x; y). If- a mapping, which is independent of the constant e, can be 
established between the vector of image values and the vector of parameters, then it will 
be possible to recover scene parameters from image data in a way that does not depend 
on illumination intensity or scene geometry. 

We now introduce the concept of an image quotient, obtained by dividing one image 
value, calculated from equation (4), by another. For a given image vector, the nth 
image quotient is given as: 

^ = — -n. m. * (5) 

Simple consideration of equation (4) shows that any quotient defined in this way will be 
invariant to a change in the parameter 8. Thus any method for the recovery of 
parameter values from image quotients will be independent of scene geometry and 
illumination intensity. 

The approach presented here is fundamentally different from most other work related to 
colour constancy. The goal in earlier work was to remove the effects of geometry and 
illmnination prior to either segmentation or pixel or object classification. Both of these 
•tasks require that pixels belonging to a particular class or object are statistically similar. 
The operations to remove the dependency on geometry and illuminating light change 
colour representation such that pixels are now clustered on the basis of the similarity of 
their intrinsic surface colouration. In contrast, the objective here is to extract 
quantitative parameters upon which the object colouration depends, not to find 
statistical similarities; Moreover, the specific filters are chosen to maximise the 
distance in the image quotient space between vectors corresponding to similar 
parameter values, as this minimises the error on the parameter value recovered from the 
colour image. 



The technique described is generally applicable to scenes in which a small number of 
parameters are required to describe all possible objects/object characteristics. In the 



formidation, the parame^^ be considered to vary continuously, Thus, the 

technique will be particularly applicable to problems where object characteristics need 
to be measured across an inoage. For example a medical imaging system may be 
required to analyse a particular tissue. The. underlying structure of the tissue will not 
vary, only specific characteristics such as thickness of the different layers (including 
zero thicloiess) or the concentration, of a particular chemical qonstituents (including zei:o 
concentration). In this situation a small parameter vector can describe all possible 
variations in the characteristics of the imaged scene! For K scene parameters the 
parameter vector is defined as: . 

K 

P = X^P/; peP (6) 

k=l 

and the space P defines all potential object characteristics. Ultimately, a mapping from 
image quotients back to the parameter vector is required, but first the forward problem 
of obtaining image quotients for* a given parameter vector is considered. A reflectance 
spectrum, corresponding to a given point within parameter space, can be described by 
the vector in M dimensional wavelength space: 

A = Am A e A . . 

where the space A defines all possible spectral reflectance functions. Hie mapping a, 
defined as 

a:P— (8) 

is introduced to denote the mapping from parameter space to wavelength space. This 
mapping gives the spectral reflectance of the object specified by the vector p. Such a 
mapping can be achieved either by a spectroscopic measurement; or by using a 
mathenaatical model which takes as input the parameters and produces a corresponding 
reflectance spectrum. Models of light propagation in different media, such as the Monte 
Carlo method or the Kubelka Munk approximation, can be used for this purpose. It 
must be possible to perform this mapping across the whole of parameter space, thus 
defining every possible spectral reflectance function. 

A digital camera with N optical filters records an JV-dimensional image vector at each 
pixel. The image vector is given as: 
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iEl (9) 

T»=l 

-wheie I describes the space of all possible image values. The .process of image 
acquisition can be considered as the projection of points in wavelength space to points 
in filter space. This projection is performed by the mapping function: 

6: A — yt. . (10) 

Equation (4) performs this mapping b in continuous form. In discrete form, the 
response of the nth optical filter, is given as: 

M 
m=l 

where Jim= Eo(X)S(/i) and the positive weights at each wavelength are given by i?^ , thus 
defining each filter response function. A digital camera effectively performs the 
mapping ib, projecting points from a large dimensional space (wavelength space) to a 
small dimensional space (filter space). With such a mapping there will be a substantial 
loss of information. However, even with this loss of information, it should be possible 
to. define the mapping in such a way that accurate information regarding the original 
parameter values can still be recovered from image data. Conditions for this mapping 
will be discussed in the following later. 

Most current image acquisition systems use an RGB system of filters. Although this 
defines a potential mapping &, it may not be the best mapping with which to recover 
parameter values from image data. However, it is known to select specific filters to 
obtain better clarity of data than that possible with an RGB system (although mainly for 
visualization or image segmentation, not for quantification). Also, in spectrometry, 
pajlicular spectral wavelengths are selected using statistical methods to improve 
quantification of components in mixtures. It will therefore be appreciated that an 
objectively defined set of optical filters is able to perform the task of recovery of 
parameters, which describe the variation in hviman skin, better than a standard RGB 
system. 



Once the vector of image values has been obtained, a vector of image quotients can be 



calculated using equation (5). The vector of image quotients is given as: 

T=f:?^ . Lei ^. . (12) 

where I describes the- space of all possible image quotients. The mapping from filter 
space to. the space of image quotients is performed after image acqmsition and will be 
referred to as mapping c, defined as: 

c:I— >L (13) 

There are many ways to define the image quotients and thus the. mapping c. For 
example, pairs of image values could be used to define image quotients as: 

■?^ = %i. n=l,2,...,^ (14) 

or a single image value could be taken as the denominator with which to calculate 
image quotients from the remaining image values, for example: 

f^ = ^ n = 2,3,...,JN?-. (15) 
^1 

At most the dimensionality, N , of the new space will be one less that that of the 
original filter space, N. This would correspond to the definition given in equation (15). 
Alternatively, if the image quotients were defined as given in equation (14), then the 
dimensionality of the new space >yill be half that of the original filter space. The aim is 
to recover a TT-dimensional parameter vector from N image quotients. Thus there must 
be at least as many image quotients as parameters, that is, N >K . 

The function/ defined -as: - , - 

/ = ao6oc . /:P — >I . . (16) 

represents the three stage mapping from parameter space to wavelength space, to image 
space, and finally to, the space of image quotients. For a given set of optical filters, it 
will be possible to perform this mapping across the whole of parameter space, provided 
that it is possible to obtain a spectrum for any given parameter vector. The inverse of 
function /is defined as: 
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and mfips jBrom the space of image quotiente directly to parametea: space. If itis possible 
to define an appropriate/"^, it will be possible to recovea: parameter values from image 
data in a way that is indepraident of illmnination intensity and scene geometry. The 
ultimate aim is to find the optimum / which maximises the accuracy of parameter 
recovery. Before a detailed discussion of this mapping is presented, it is important to 
emphasise that the fonn of the function / will depend on the mappings a, b and c. 
Although mapping a is fixed for a given problem, mapping b will vary with the choice 
of optical filters and mapping c will vary depending on how the image quotients are 
■ defined. 

Any mapping function which is to map from the space of image quotients (I -space) to 
parameter space (P-space) must be 1 to 1. That is, for a given point in P-space, th^e 
must be a corresponding unique point in I -space and vice-versa. If this is not the case, 
ambiguity will arise as it could be possible to recover more that one set of parameter 
values from a given vector of image quotients. Once this condition has been 
established, it is necessary to consider the error associated with parameter recovery as, 
using a digital camera, it will only be possible to obtain image values to within a given 
uncertainty. This will introduce an uncertainty into the final recovery of the parameter 
vector. There could also be an error associated with the prediction/measurement of a 
spedtrum from the parameter vector. For simplicity the analysis presented here will be 
restricted to problems in which the error associated with the spectral measurement can 
be neglected. 

Initially, the problem where one parameter is sufficient to describe all variation in an 
imaged scene will be analysed. The metiiodology will then be extended to problems 
where the number of parameters is greater than one. 

Consider the case where one image quotient (two image values) is used to recover a 
single parameter value. Figure 1 illustrates a function /which gives the image quotient 
as a function of the parameter p. It is clear that in order to satisfy the 1 to 1 condition, 
the curve must not have any turning points: that is, it must increase or decrease 



monotonically in the appropriate range of p. Mathematically this is expressed as: • , 
>0 Vp€P. - ' - (18) 



dp 



Measurement of an image quotient value ig, corresponding to a parameter value po, is 
now considered. Associated with acquisition of each image value is an uncertainty due 
to camera error. It is straightfprward to show, using standard error analysis, that the 
error associated with an image quotient T, which has been calculated from the two 
image values i] and 13, is given as: 

where A£ is the camera uncertainty. This error has been shown on iJie ordinate of the 
graph in Figure 1. If the derivate of /is non-zero in some neighbourhood of po then it is 
possible to approximate this function linearly. Assuming the error Ai to lie within this 
neighbourhood, the corresponding error in the parameter value is given as: 

■Ap^Ai/^. - • • (20) 

Thus, it is possible -to obtain a value for the error A/?, associated with parameter 
recovery, at any point in P-space. An optimisation criterion can then be defined based 
on some measure of this error. For most applications it will be necessary to minimise 
the .error equally across the whole of P-space. For others it may be that high accuracy 
parameter recovery is required within a certain range of parameter values. For example, 
in a medical image application, imaged tissue could be deemed pathological once a 
characterising parameter changes beyond a threshold level. This would need to be 
accounted for with" some form of weighting in the optimisation criterion. It is 
interesting to note that in order to minimse A/?, it is necessary to maximise the 
magnitude of the derivative given in equation (18). This wiU ensure that any search, 
carried out to find an optimum /, will tend to move towards regions of search space 
where the 1 to 1 condition is satisfied. 



In theory it is possible to recover the parameter using more than one image quotient. In 
this case it wiU be necessary to calculate the error associated with parameter recovery 
for each of the image quotients and select the one, at each point in P-space, which has 
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the smallest associated error (Ap). It may be that the optimisation procedure gives a 
single image.quotient which perfomis better than any other across the whole of P-space. 
In this situation there is" no benefit to using more Biat one image quotient. 

The analysis is now extended to the general problem where the recovery of a Z*- 
dimensional parameter vector is required from an N dimensional vector of image^ 
quotients. Initially the analysis will be restricted to the case where N - K and will then 
be extended to include situations where N >K. As discussed earlier, if N <K, dien it 
is not possible to recover Z-dimensional data from an N dimensional measurement 



The mapping function f , defined as: 

i=f(p) 

* must now be considered a vector valued function of a vector variable. In the following 
analysis specific results from differential geometry will be used. For further details the 
reader is directed to for example M. M. Lipschutz, Differential geometry (McGraw-Hill 
Book Company, New York, 1969). To establish whether the function/provides a 1 to 1 
relationship, it is first necessary to consider the behaviour of the determinant of the 
Jacobian matrix, simply referred to as the Jacobian. This is defined as: 



dpi dp2 



dpi dp2 



dh 

dPk 
dPh 



(22) 



The Jacobian can be considered the multidimensional equivalent of the one dimensional 
derivative given in equation (18). The inverse function theorem states ttiat, if the 
Jacobian is non-zero at a point po in P-space, then there exists a neighbourhood around 
po where the function f can be approximated linearly as 

f (p) = f (Po) + ftf (Po)(p - Po) 



where is the differential of f and is given as: 



(23) 



(24) 



It follows that in this neighboudiood the function f provides a 1 to 1 relationship. Thus, 
if it is possible to establish that the Jacobian is strictly positive or strictly negative 
throughout the whole of P^space, the function f will be 1 to 1 eyerywhete. Once this 
. condition has been established, it is necessary to consider how the error associated with 
image acquisition maps under f ^ to give the corresponding error in parameter recovery. 
The error associated with each image quotient is calculated using equation (19). The 
combination of errors maps out a hypervolume in I -space, centred on the point io. This 
has been illustrated in Figure 2 for the case of a 2D P-space, where an ellipse is 
obtained (or a circle if the errors are equal). An ellipsoid is obtained in 3D space and a 
hyperellipse in higher dimensions. Although the following analysis will be based on a 
2D P-space, the arguments are equally valid in higher dimensions. 

The ellipse in I -space represents all possible image quotient vectors which could 
correspond to a camera measurement i = Iq . It is assumed that the region of error lies 
within the neighbourhood of i = ij, where the mapping function f can be approximated 

linearly. Thus, under the mapping f ^ the ellipse in I -space maps directly to another 
ellipse in P-space. This new ellipse is centred on the point p = po and represents all 
possible parameter vectors which could be recovered from the vector of iniage quotients 
i = io • The error associated with parameter recovery is obtained by considering the 
worst case scenario: that is the point within the ellipse in P-space which is at the 
maximum distance from the point p = po. This maximum distance must- be calculated 
separately for each component, pk, of the parameter vector to obtain the error associated 
with recovery pf each individual component. To calculate these errors it is necessary to 
consider how the ellipse is transformed under the mapping f \ which is linear provided 
the Jacobian is non-zero. 

> Under a Unear mapping the ellipse will be translated, scaled and rotated. The translation 
associated with the linear mapping defines ttie point p = po which is mapped to from the 
point i = to . The two other transformations, scaling and rotation, are best understood by 

considering how a vector di = <ff in I -space, maps under to give a corresponding 
vector dg in P-space. ■ The vector dp can be calculated, from the inverse form of 
•equation (24) which, in matrix form, is given as: 
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dp==J-^cfi (25) 

where J denotes the.Jacobian matrix. Note that JT^ exists only if the Jacobian in non- 
zero. This must be the case if the 1 to 1 condition is to be satisfied. 



The vectors A and B correspond to the major and minor axes of the ellipse in I -space 
and are given as: 



/ _ \ 




Ail 








H 











(26) 









Apf 


A' = 


■ B' = 











Under the mapping f ^ these vectors map to the vectors A' and B' which correspond to 
the major and nrinor axes of the ellipse in parameter space. Solving equation (25) for 
each of these vectors gives: 



(27) 



where Apj"* and tq>^ are the components of the vector A' in the direction of p\ and 
respectively. Similarly, Apf and Ap^ are the components of the vector B' in the 
direction of p\ and pn respectively. To calculate the error in each component of the 
parameter vector it is necessary to consider the worst case scenario. It can be seen from 
Figure 2 that this corresponds to taking the maximxmi of Lp^ and Ap^ as the error in pi 
and taking the maximum Ap^ and Kpl as the error in po. This error can be specified 
by a vector Ap and can be calculated for any given point in parameter space. With this 
measure of the accuracy of parameter recovery across the whole of parameter space, it 
is possible to define an optinodsation criterion. This could simply be based on a sum of 
the errors at every point in P-space or could be chosen to favour accuracy of recovery of 
a' subset of the original parameters. Once this optimisation criterion has been defined, a 
search can be used to find the optimum mapping function f . It is important to note that, 
although the above discussion is based on a 2D parameter space, the methodology is 
equaUy applicable to any iT-dimensional parameter space. 



An algorithm for the implementation of the proposed methodology is given as follows: 



. . i . Establish a suitable search space from a parameterisation of mappings b and c. 
2. For a given mapping function f calculate the vector of image quotients for 
each point within a discretised parameter space. " ■ 

" 3. For each point, check that the Jacobian is either strictly positive or strictly 
negative across the whole of parameter space. If this condition is held then 
compute, the inverse of the Jaco]bian matrix. If not then jetum to step 1 and 
define a new mapping function f . 

4. Using equation (19), calculate the error associated with image acquisition for 
each vector of image quotients. 

5. From the image quotient vector error calculate the maximum possible error in 
each component, pj^, of the parameter vector across the whole of parameter 
space. 

6. Use the vector of parameter errors at each point within parameter space to 
measure the accuracy of parameter recovery. 

7. Repeat steps 2-6 with some optimisation technique which enables an optimum 
mappiiig function f to be determined. 

It is fairly straightforward to extend this, methodology to the case in which N > K : that 
is, where there are more image quotients than parameter values. Initially every possible 
K-dimensional subspace of image quotients will need to be defined from the original- 
N -dimensional space of image quotients. It wiU then be necessary to go through the 
above procedure for each potential subspace and obtain the vector of parameter errors at 
each point within parameter space. To achieve the maximum' possible accuracy the best 
Ap must be selected at every location within parameter space. Thus every point in P- 
space will be linked to a specific image quotient combination. It will then be necessary 
to link every region of the original N -dimensional space of image quotients to the 
particular subspace of image quotients which should be used for parameter recovery. It 
is important to note that it is necessary to recover the whole parameter vector at each- 
point within a particular i?-dimensional subspace of image quotients. It is not 
possible to attempt to improve the accuracy of the system by recovering different 
components of the parameter vector from different jfiT-dimensional subspaces of image 
quotients. This is mathematically invalid. 

The mapping function f is a composite function of three separate mappings. Although 
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the first mapping a, from parameter space to wavelength space, is fixed for a given 
problem, mappings b and c can vary depending on the choice of optical filters and 
definition of image quotients. Thus, to define an appropriate search -space it is 
necessary to parameterise -mappings b and c. Mapping fe, which represents image 
acquisition, is defined by the positive NxM matrix i?^, given in equation (11). 
Typically this matrix will contain many elements and an appropriate parameterisation 
should be based on typical filter response functions. For example, the position of the 
central wavelength and a measure of width could be used to define a Gaussian shape. 

Parameterisation of the mapping function c wiU be fairly straightforward as there are 
only a liroited number of ways of combining image values to produce independent 
image quotients. In .some applications the form of this mapping may be fixed apriori. 
Thus, it will not increase the overall dimensionality of the search space. 

An optimisation method should search the whole space of possible mappings using the 
optimisation criterion outUned in the previous section. One technique which is ideally 
suited to this type of search is a genetic algorithm, GA, {see T. Back and H. P. 
Schwefel, "An overview of evolutionary algorithms for parameter optimisation," 
Evolutionary Computation 1, 1-23 (1993)} as it is straightforward to define a fitness 
function which measures the accuracy of parameter recovery. Genetic algorithms have 
been shown to work well on a wide range of problems with objective functions that do 
not possess "nice" properties such as continuity, differentiability or satisfaction of the 
Lipschitz Condition {see L. Davis, The handbook of genetic algorithms (Van Nostrand 
Reingold, New York, 1991), and D. Goldberg, Genetic algorithms in search, 
optimization and machine learning (Addison-Wesley, London, 1989)}. 

The above techniques will now be further exemplified by considering their application 
to the analysis of normal skin composition. Firsfly, the prediction of spectral 
reflectance is considered. 

In order to perform mapping a it is necessary to have either a mathematical model 
which can predict spectral reflectance for a given set of parameter values or some 
technique for measurement of the appropriate spectrum. For this application we use the 



mathematical^ i^^ developed by Cotton and Claridge {see S. D. Cotton and E 
•Claridge, ^Developing a predictive model' of human skin colouring," Proc, of SPIE 
Med. Lnag. 2708, 814-825 (1996)}. With this modehit is' possible to predict tiie 
spectral reflectance for a given set'of parameters. An outiine of the model is now given. 

Skin can be considered to be the four-layer structure depicted in Figure 3. A negligible 
amount of light is reflected from the surface of the skin, thus the surface term in 
equation (2) can be neglected. Although not absorbing any radiation, the stratum 
comeum scatters the incoming light in all directions. light which penetrates this layer 
can thus be considered diffuse. In the epidermis light is absorbed by the pigment 
melanin. The absorption at each wavelength can be calculated using the- Lambert-Beer 
law and will depend on the product of the melanin extinction coefficient and the 
pigment concentration. After passing through the epideraiis the light is^both scattered 
and absorbed by the papillary dermis. The absorption results from the presence of 
-blood and scattering from the underlying collagen structure. The simple Kubelka- 
Munk light tiieory {P. Kubelka and F Munk, "Bin Beitrag zur Optik der Farbanstriche", 
"Z. Tech. Opt" 11, 593-611 (1931)} can be used to model the interaction of Ught with 
the papillary demiis as the necessary condition of diffuse incident illumination is 
satisfied. Any light which passes through the papillary dermis into the recticular dermis * 
can be neglected as no significant backscattering occurs in this layer. Using this two- 
layer light transport model it is possible to obtain the remitted spectra for given 
concentrations of melanin and blood. A more detailed description of 'this model can be 
found in S. D. Cotton and E Claridge, "Developing a predictive model of human skin 
colouring", Proc. of SPIE Med. Imag. 2708, 814-825 (1996). 

For a given papillary dermal thickness, changes in melanin and blood characterise all 
histological variation and thus define a 2-D parameter space for healthy skin. To-cairy 
out the optimisation procedure described above it is necessary to discretise parameter 
space. This is done at equal intervals to define 10x10 points, each of which 
corresponds to a spectrum generated by the mathematical model. For simplicity, 
concentration values wiU be denoted by a number between 1 and 10. Figures 4a and 4b 
show how the remitted spectrum changes as melanin and blood are varied respectively. 
•With a change in melanin concentration; the intensity of the whole spectrum is seen to 
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decrease, with a more pronounced change in the blue region. As the blood 

concentration is decreased the most significant reduction in intensity is observed in the 

green region, the resulting shape reflecting "the two absorption maxima of 

oxyhaemoglobin, a blood bom pigment " . "* 

« 

To define a siiitable search space it is necessary to parameterise the mappings b and c. 
A'parameterised form of b is chosen to define a typical interference filter. This is 
modelled as a square profile with Gaussian decay at each side. Two parameters are 
required to specify this shape: the central wavelength and a full width half maximum 
(FWHM). Optimisation is carried out for three such filters, defining a 6-D search space. 
With three filters giving three.image values, h; 12 and 13, the only possible definition of 
image quotients, if we assume ii = 13 is equivalent to 13 = iu is given as: 

= (28) 
In this instance the mapping c does not increase the dimensionality of the search space. 

The optimisation procedure was implemented following the algorithm given above. 
Initially the vector of image quotients was_ calculated for every point within the 
. discretised parameter space. This was done using the mathematical model to perform 
mapping a, the parameterised form of matrix to perform mapping b and the 
equations (28) to perform mapping c. The derivative of each image quotient, with 
respect to each parameter, was obtained at each point within discretised parameter space 
using three-point finite difference approximations. The Jacobian matrix was then 
constracted at every point within parameter space, and providing its detemiinant was 
* non-zero everywhere, the inverse calculated. If this condition was violated then a new 
mapping f was defined. The errors associated with image acquisition were then 
calculated using equation (19). The absolute value of the error in each image value will 
vary depending on the camera gain setting. Although this constant will npt affect the 
mapping f , it must be estimated in order to calculate the effective camera- error. For this 
appUcation it was taken to be 0.78% of the maximxmi value of all the image values 
across parameter space. This corresponds to a camera which has been set to give a 
' maximum reading for the largest spectral reflectance and a camera error of two grey 
scale levels in an 8-bit representation. 



Using the procedure outlined above the error associated with parameter recovery in both 
melanin and blood was obtained for each point within the discretised parameter space. 
In order to "find an optimum f, it is necessary to minimise the errors in recovery of both" 
melanin and blood across the whole of parameter space. Thus the fitness function for 
the GA was taken to be the sum of ttie errors in both melanin and blood. This procedure 
was implemented in'matlab*™ using a standard GA to search the space of available 
mappings. 

The boundaries of the search space were chosen such that the central wavelength was 
constrained to lie in the visible region (400nm-700imi) and such that the widths of the 
filters were allowed to vary from a FWHM of 25 to 200nm. Although it is now possible 
to engineer almost any shape of interference filter, this corresponds to an economically 
viable range of such filters. 

Although it was originally assumed that an image quotient defined as ii I z's would be 
equiv^ent 'to h I ii, the results of the GA search showed that this was not the case. The 
search was intitialised for a random seed and, although the same central wavelengths 
were always obtained, different filters were selected corresponding to 13 defined in 
equation (28). Further investigation showed that these local maxima in the search space 
corresponded to differing distributions of errors both, across parameter space and 
between the two parameters.* This is because the fitness function, or measure of 
accuracy, was defined as the sum the errors across parameter space for both melanin 
and blood. Thus, a loss of accuracy in one parameter could be compensated for with an 
increase. in the other. It may be that, with a more exact specification of the error 
distribution in the fitness function, it woxxld be possible to obtain the same results for 
every GA search. . .. 

Figure 5 shows a filter combination which gave a similar error in the recovery of both 
melanin and blood. The image quotients were calculated by dividing the filter centred 
at A. = 473nm and X = 560nm by the response of the filter centred at A. = 700nm. To 
understand why these specific filters were selected it is necessary to analyse the spectral 
curves shown in Figure 4. The filters centred at X = 473nm and X = 560nimmi 
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coirespond to spectral locations where there is a large change in intensity with the 
parameters melanin and blood respectively; A third filter was then required in a region 
of the spectrum in which the remitted light which was either significantly less or 
significantly more than that of the other two filters. The filter centred at X = 700mri was 

• chosen as it always gave the largest response at any point within parameter space. This 
ensured that the derivatives of each image quotient decreased monotonically across the 
whole of parameter space. The Jacobian, calculated firom these derivatives, was strictly 
positive across the whole of parameter space. It is interesting to note that some 
alternative filter combinations gave 'jacobians which were strictly negative across 
parameter space, corresponding to alternative local maxima in the search sp£w». If two 
filters are chosen, to define an image quotient, which vary similarly across parameter 
space, there will be minimal change in that image quotient and thus it wiU be of limited 
value for parameter measurement 

Table 1 and Table 2 give the errors associated with recovery of melanin and blood 
respectively at every point vyithin the discretised parameter space, for the filter 
combination given in Figure 5. These errors have been calculated assuming a camera 
uncertainty of two grey scale levels in an 8-bit representation that is, an error of (2/255). 
The errors can be seen to increase as both the melanin and blood concentrations 
increase. This is due to the decrease in the magnitude of the final remitted spectrum as 
the pigment concentrations increase. Although this absolute error increases, the 
fractional error at each point acmally decreases and it is nearly always possible to 

• differentiate points in parameter space as the error rarely rises above 0.5. 

it has been demonstrated that, using an objectively defined set of optical filters, it is 
possible to recover scene parameters from image data in a way which is insensitive to 
geometry and incident illumination. In the example problem, discussed above, the error 
associated with this parameter recovery was found to be relatively small. The 
invariance of this mapping means that the technique will be particularly applicable to 
medical image applications where there is significant curvature of the surface of the 
imaged tissue, such as near a joint. It will also be unnecessary to calibrate the camera to 
determine the intensity of the incident Ught. This could help to significantly increase 
the speed of image acquisition and later processing. 



The methodology set out here has been developed for- a measurement task, where the 
scene parameters are known to vary continuously. The technique can be also be applied 
to problems of "recognition, where it is necessary to differentiate discrete objects based 
on some measure of their spectral reflectance. This approach has been discussed in the 
article G. Healey, "Using .colour for geometry-insensitive segmentation," J. .Opt. Soc. 
Am. A 6, 920-937 (1989) who used the idea of normalised colour to identify different 
regions of normalised colour space • corresponding to different metal and di-electric 
materials. This enabled geometry-insensitive segmentation of an imaged object 
comprised of a number of different materials. 

It will be appreciated that in order to implement the proposed methodology, a look-up 
table should be established between all possible image quotients arid scene parameters. 
Although this may be time consuming, it is only necessary to carry out this procedure 
once. Once established, this look-up table will ensure no significant processing after 
image acquisition, making this technique particularly suitable to real-time applications. 

Figure 6 is a flow diagram showing the key steps in the method described above. 

It will be appreciated by the person of skill in the art that various modifications may be 
made to the aLbove described embodiments without departing from the scope of the 
present invention. " • 
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Claims: 

1. 'A method of determining a property or properties of each of a set of filfeii, 
which filters are used to select specific wavelength ranges in" a system which relies upon 
a speclxal analysis of remitted, emitted, and/or transmitted Ught to quantify a parameter 
or parameters of an object or material, the method comprising the steps of: 

1) defining a set of potential filter properties 

2) defining one or inore image quotients, the or each image quotient for a region 
being, obtained by dividing the quantified output of a given filler for that region, the 
"image value" for that filter, by another image value for that same region; 

3) for an object or material to be analysed and for said defined set of potential 
filter properties and for said, given image quotients, olptaining a function mapping points 
in parameter space to points in image quotient space; 

4) determining whether the mapping function provides a 1:1 correspondence 
between points in parameter space and points in image quotient space; and 

■ 5) if the mapping function does not provide a 1 : 1 correspondence, rejecting said 
potential filter properties, repeating steps 1) to 4) and, if the mapping function does 
provide a 1:1 correspondence accepting the potential filter properties as a candidate set 
of filter properties. 

2. A method according to claim 1 and comprising, for a set of potential filter 
properties accepted as a candidate set of filter properties, detennining the accuracy of 
parameter recovery obtained using said mapping function and determining whether or 
not the accuracy is sufficient or matches some other criterion. 

3. A method according to claim 2, wherein if the accuracy is sufficient or matches 
said criterion, the candidate filter properties are adopted and if the accuracy is not 
sufficient or does not match said criterion, steps 1) to 5) are repeated for a different set 
of filter properties. 

4. A method according to claim 2 or 3 , wherein the accuracy is determined by: 
a) calculating the error associated with image acquisition for each vector 

of each image quotients; 




. b) from the image quotient vectpr error, calculating the ma?dmum possible error 
in each component of the parameter vector across the whole of parameter space; 
an d . - . , - 

.c) Using the vector of parameter errors at each point Within parameter space to " 
measure the. accuracy of parameter recovery. 

5. A method according to claim 1 or 2 and comprising repeating steps 1) to 5) for a 
multiplicity of sets of potential filter properties to identify a plurality of candidate filter 
property sets, detemuning for each candidate set an error value representing the 
accuracy of parameter recovery obtained using the corresponding mapping function, 
and using said candidate set as a basis for determining a preferred set of filter properties. 

6. A method according to claim* 5 and comprising using a genetic algorithm to 
determine a preferred set of filter diaracteristics using said candidate set. 

7. A method according to' any one of claims 1 to 3 and comprising using a gradient 
descent algorithm to select an optimal set of filter characteristics, the starting point for 
the algorithm being a first candidate set of filter properties identified in step 4). 

8. A method according to any one of the preceding claims, wherein for each filter 
the method of the present invention is used to determine the centre wavdiength of the 
filter. 

9. A method according to any one of the preceding claims, wherein the method is 
used to determine the fiill width half maximum OPWHM) of the filter characteristic. 

10. A method according to any one of the preceding claims, wherein step -3) 
comprises constructing a Jacobian matrix for the mapping function with respect to said 
parameter(s), and obtaining the determinant of that matrix. 

11. A method of analysing body tissue or fluid and comprising using the method of 
any one of the preceding claims. 
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12. A method of determining a' property or jwoperties of each of a set of filths 
substantially as hereinbefore described with reference to the accompanying drawings. 

13. Apparatus for analysing an object or material having means for conducting a 
spectral analysis of remitted, emitted, and/or transmitted Hght to quantify a parameter or 
parameters of an object or material, "the apparatus comprising a plurality of filters for' 
spMtting said light into respective components, the filters having properties obtained by 
using the method of the above first aspect of the invention. 



ABSI^ 

Method and Apparatus for Quantifying INfoterial or ObjectProperties 

A method of determining a properties of each of a set of filters, whieh filters are used to 
select specific wavelength ranges in a system which relies upon a spectral analysis of 
remitted light to quantify parameters of an object or material. The method comprises 
defining one or more image quotients. For an object or material to be analysed and for a 
given set of potential filter properties, a fimction mapping points in parameter space to 
points in image quotient space is obtained It is then determined whether the mapping 
function provides a 1:1 correspondence between points in parameter space and points in 
image quotient space. If the mapping function does not provide a 1:1 correspondence, 
said potential filter properties are rejected and, steps 2) and 3) repeated If the mapping 
function does provide a 1:1 correspondence, the potential filter properties are accepted 
as a candidate set of filter properties. 
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Figure 4b 
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Define discrete representation of spaces P, A, I and I 
and dejBner a mapping function a: P A 



Define filter properties (a mapping function b: A — > I) and 
define image quotients (a mapping function c: I — > I) 
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For the given definition of filter properties and the given 
definition of image quotients obtain a mapping function 

f:P-^ I 



Using the given mapping function f calculate the vector of 
image quotients for each point (parameter vector) p within 

the discretised parameter space P 



For each point p in P compute the Jacobian df/dp 



Amend filter properties 
and/or image quotients 




NO 



YES 



For each vector of image quotients calculate the error 
associated with image acquisition 



From the image quotient error calculate the maximum error in 
each component pk.of the parameter vector p_ across the whole 
. of parameter space P 



Measure the accuracy of parameter recovery using the vector 
of parameter errors at each point p in parameter space P 



Amend filter properties 
and/or image quotients 



Accuracy of 
parameter recovery optimal and/or 
^ sufficient? 

NO 

-< YES 



Figure 6 



This Page is Inserted by IFW Indexing and Scanning 
Operations and is not part of the Official Record. 

BEST AVAILABLE IMAGES 

Defective images within this document are accurate representations of the original 
documents submitted by the applicant. 

Defects in the images include but are not limited to the items checked: 

□ BLACK BORDERS 

□ IMAGE CUT OFF AT TOP, BOTTOM OR SffiES 

□ FADED TEXT OR DRAWING 

□ BLURRED OR ILLEGIBLE TEXT OR DRAWING 

□ SKEWED/SLANTED IMAGES 

^ COLOR OR BLACK AND WHITE PHOTOGRAPHS 

□ GRAY SCALE DOCUMENTS 

Q^ilNES OR MARKS ON ORIGINAL DOCUMENT 

□ REFERENCE(S) OR EXHIBIT(S) SUBMITTED ARE POOR QUALITY 

□ OTHER: 

IMAGES ARE BEST AVAILABLE COPY. 
As rescanning these documents will not correct the image 
problems checked, please do not report these problems to 
the IFW Image Problem Mailbox. 



